Confinement-Shear Lattice Model for Concrete Damage in Tension and Compression: II. Computation and Validation

The concrete material model developed in the preceding Part I of this study is formulated numerically. The new model is then verified by comparisons with experimental data for compressive and tensile uniaxial tests, biaxial tests, and triaxial tests, as well as notched tests of mode I fracture and size effect. DOI: 10.1061/~ASCE!0733-9399~2003!129:12~1449! CE Database subject headings: Concrete; Microstructure; Fractures; Damage; Softening; Failures; Materials tests; Computer analysis; Nonlinear analysis; Particle interactions; Particle distribution. of the element is Introduction In Part I of this study, a new three-dimensional constitutive model for concrete has been formulated as a three-dimensional lattice connecting the centers of randomly distributed aggregate particles. In the present Part II, a numerical algorithm for this model will be presented and the model will be validated by comparing numerical results and experimental data found in the literature for typical uniaxial, biaxial, and triaxial tests, as well as notched tests of model I fracture and size effect. Only short-time loading, for which creep is unimportant, and loading rates low enough for dynamic effects to be absent, will be considered. All the definitions and notations from part I will be retained. Numerical Implementation and Stability Conditions Concrete is here idealized as an assembly of particles whose interaction can be described by a two-node finite element ~Cusatis et al. 2003!. Like the finite element method in general, this leads to the following matrix equation of motion for the structure: MQ̈1CQ̇1P~Q!5F~ t ! (1) Here M5mass matrix of components M i ; P5vector of internal forces of components Pi , which are obtained by assembling the contributions from all the finite elements ~all struts!; Q5the vector of kinematic variables Qi ~displacements and rotations of all Graduate Student, Dept. of Structural Engineering, Technical Univ. ~Politecnico! of Milan, Milan 20133, Italy. E-mail: [email protected] McCormick School Professor and W. P. Murphy Professor of Civil Engineering and Materials Science, Northwestern Univ., Evanston, IL 60208. E-mail: [email protected] Professor of Structural Engineering, Dept. of Structural Engineering, Technical Univ. ~Politecnico! of Milan, Milan 20133, Italy. E-mail: [email protected] Note. Associate Editor: Stein Sture. Discussion open until May 1, 2004. Separate discussions must be submitted for individual papers. To extend the closing date by one month, a written request must be filed with the ASCE Managing Editor. The manuscript for this paper was submitted for review and possible publication on August 26, 2002; approved on February 21, 2003. This paper is part of the Journal of Engineering Mechanics, Vol. 129, No. 12, December 1, 2003. ©ASCE, ISSN 07339399/2003/12-1449–1458/$18.00. JOURN aggregates!; F(t)5given load history; and subscripts i51,2,...,N number the generalized displacement components. The damping matrix C, of components Ci , is introduced only for the sake of generality and is assumed zero in computations. Eq. ~1! can be integrated explicitly or implicitly. The implicit methods allow much larger time steps but require considerably more computational effort per time step than the explicit methods. The explicit methods have severe limitations on the time step length because of the condition of numerical stability, but have the advantage of much smaller memory requirements. In the present three-dimensional numerical simulations, there are thousands of unknowns even for the smallest test specimen, and so computer memory limitations and running time are important. Therefore, the explicit method, based on a central difference algorithm, is selected ~it allowed running all the present problems on a personal computer with a 750 MHz Pentium III processor!. The central difference approximation of ~1! is Q̇i 5aQ̇i 1b~Pi 2Fi !/M i (2) for unconstrained degrees of freedom and Q̇i 5 ḋ i 5 ḋ i~ t n11/2!, Ri 5Pi 1M i~ ḋ i 2a ḋ i !/b (3) for constrained degrees of freedom ~see the Appendix!; ḋ i(t)5prescribed velocity history for the ith degree of freedom; Ri 5work-conjugate reaction of the constraint; n51,2,...,M are the labels for discrete times tn subdividing time t into intervals Dt . The parameter values are: a50, b5Dt/2 for n50, and a 5(12v)/(11v), b5Dt/(11v) for n.0, where v5CDt/2, Dt5time step, C5damping coefficient, and C5CM ~which was considered as zero in all calculations!. The numerical stability condition simply is Dt,2/v ~Belytschko et al. 2000!, where v represents the highest natural frequency of the system ~Bažant and Cedolin 1991!. Because the longitudinal displacement field throughout the element is linear, the frequency of free vibrations in the longitudinal direction is vN 2 54EN /(rl ) ~Belytschko et al. 2000!. For the shear free modes, a similar formula can be derived. In two dimensions ~the extension to three dimensions being straightforward! the kinetic energy associated with the shear deformation AL OF ENGINEERING MECHANICS © ASCE / DECEMBER 2003 / 1449 KT5 1 2 E0 l1 rA@ v̇~x1!#dx11 1 2 E0 l2 rA@ v̇~x2!#dx2 (4) where v(x1) and v(x2)5transversal velocities of particles 1 and 2, respectively; v(x1)5v11x1q1 , v(x2)5v21x2q2 , and v̇(x1)5 v̇11x1q̇1 , v̇(x2)5 v̇21x2q̇2 ~Cusatis et al. 2003!. By substituting these expressions into Eq. ~4!, we get KT 5q̇MT c q̇/2, where q̇5$v̇1 ,q̇1 , v̇2 ,q̇1% and MT c 5rAF l1 l1/2 0 0 l1/2 l1/3 0 0 0 0 l2 2l2/2 0 0 2l2 /2 l2 /3 G (5) The lumped mass matrix can be obtained by simply neglecting the terms out of the diagonal, as usually done for beam elements ~Belytschko et al. 2000!; MT5rA diag@l1 ,l1 /3,l2 ,l2 /3# . The elastic energy associated with these free modes can be expressed as WT5*0 AET«Tdx/2. Since the shear strain is «T5(v22v1 2l2q22l1q1)/l ~Cusatis et al. 2003! one gets WT5qKTq/2, where q5$v1 ,q1 ,v2 ,q1% and